Consider a coin-tossing experiment with outcomes $x \in\{0,1\}$ (tail and head, respectively) and let $0\leq \mu \leq 1$ represent the probability of heads. The data generating distribution for this model can written as a Bernoulli distribution:
$$ p(x|\mu) = \mu^{x}(1-\mu)^{1-x}$$
Note that the variable $x$ acts as a (binary) selector for the tail or head probabilities. Think of this as an 'if'-statement in programming.
We recognize the $(\alpha_k)$'s as prior pseudo-counts and the Dirichlet distribution shows to be a conjugate prior to the categorical/multinomial:
$$\begin{align*} \underbrace{\text{Dirichlet}}_{\text{posterior}} &\propto \underbrace{\text{categorical}}_{\text{likelihood}} \cdot \underbrace{\text{Dirichlet}}_{\text{prior}} \end{align*}$$
A related distribution is the distribution over count observations $D_m=\{m_1,\ldots,m_K\}$, which is called the multinomial distribution,
$$p(D_m|\mu) =\frac{N!}{m_1! m_2!\ldots m_K!} \,\prod_k \mu_k^{m_k}\,.$$
Now let's proceed with learning the parameters for a model for $N$ independent-and-identically-distributed (IID) rolls of a $K$-sided die, based on observed data set $D=\{x_1,\ldots,x_N\}$.
As above, consider the following model assumptions for $N$ tosses with a $K$-sided die with parameters $\mu = (\mu_1,\mu_2, \ldots,\mu_K)$.
$$\begin{align} p(D|\mu) &= \prod_{n=1}^N \mathrm{Cat}(x_n|\mu) = \prod_{k=1}^{K} \mu_k^{m_k} \tag{likelihood}\\ p(\mu|\alpha) &= \mathrm{Dir}(\mu|\alpha) = \frac{1}{B(\alpha)} \prod_{k=1}^{K} \mu_k^{\alpha_k -1} \tag{prior} \end{align}$$
where $B(\alpha) = \frac{\prod_k \Gamma(\alpha_k)}{\Gamma(\sum_k \alpha_k)}$ is known as the Beta function.
Work out both the model evidence and the posterior distribution for $\mu$.
🎯 Key concept
Discrete event outcomes are typically represented via one-hot encoding, in which each outcome corresponds to a unique binary indicator vector.
Let's apply what we have learned about the loaded die to compute the probability that we throw the $k$-th face at the next toss.
$$\begin{align*} p(x_{\bullet,k}=1|D) &= \int p(x_{\bullet,k}=1|\mu)\,p(\mu|D) \,\mathrm{d}\mu \\ &= \int_0^1 \mu_k \times \mathcal{Dir}(\mu|\,\alpha+m) \,\mathrm{d}\mu \\ &= \mathrm{E}\left[ \mu_k | D\right] \\ &= \frac{m_k + \alpha_k }{ N+ \sum_k \alpha_k} \end{align*}$$
(You can find the mean of the Dirichlet distribution $\mathrm{E}\left[ \mu_k \right]$ at its Wikipedia site).
This result is simply a generalization of Laplace's rule of succession.
The posterior for $\{\mu_k\}$ can be obtained through Bayes rule:
$$\begin{align*} p(\mu|D,\alpha) &\propto p(D|\mu) \cdot p(\mu|\alpha) \\ &\propto \prod_k \mu_k^{m_k} \cdot \prod_k \mu_k^{\alpha_k-1} \\ &= \prod_k \mu_k^{\alpha_k + m_k -1}\\ &\propto \mathrm{Dir}\left(\mu\,|\,\alpha + m \right) \tag{B-2.41} \\ &= \frac{\Gamma\left(\sum_k (\alpha_k + m_k) \right)}{\Gamma(\alpha_1+m_1) \Gamma(\alpha_2+m_2) \cdots \Gamma(\alpha_K + m_K)} \prod_{k=1}^K \mu_k^{\alpha_k + m_k -1} \end{align*}$$
where $m = (m_1,m_2,\ldots,m_K)^T$ is the count vector.
Setting the derivative of $\tilde{\mathrm{L}}(\mu)$ to zero yields the sample proportion for $\mu_k$
$$\begin{equation*} \nabla_{\mu_k} \tilde{\mathrm{L}}(\mu) = \frac{m_k } {\hat\mu_k } - \lambda \overset{!}{=} 0 \; \Rightarrow \; \hat\mu_k = \frac{m_k }{N} \end{equation*}$$
where we get $\lambda$ from the constraint
$$\begin{equation*} \sum_k \hat \mu_k = \sum_k \frac{m_k} {\lambda} = \frac{N}{\lambda} \overset{!}{=} 1 \end{equation*}$$
(We insert this slide only to alert you to the difference between using one-hot encoded outcomes $D=\{x_1,x_2,\ldots,x_N\}$ as the data, versus using counts $D_m = \{m_1,m_2,\ldots,m_K\}$ as the data. When used as a likelihood function for $\mu$, it makes no difference whether you use $p(D|\mu)$ or $p(D_m|\mu)$.)
Of course, we shouldn't have to go through the full Bayesian framework to get the maximum likelihood estimate. Alternatively, we can find the maximum likelihood (ML) solution directly by optimizing the (constrained) log-likelihood.
The log-likelihood for the multinomial distribution is given by
$$\begin{align*} \mathrm{L}(\mu) &\triangleq \log p(D_m|\mu) \propto \log \prod_k \mu_k^{m_k} = \sum_k m_k \log \mu_k \end{align*}$$
$$\begin{align} \overbrace{\prod_{k=1}^{K} \mu_k^{m_k}}^{\text{likelihood }p(D|\mu)} \cdot \overbrace{\frac{1}{B(\alpha)} \prod_{k=1}^{K} \mu_k^{\alpha_k -1}}^{\text{prior }p(\mu|\alpha)} &= \frac{1}{B(\alpha)} \prod_{k=1}^{K} \mu_k^{m_k + \alpha_k -1} \\ &= \frac{B(m+\alpha)}{B(\alpha)} \frac{1}{B(m+\alpha)}\prod_{k=1}^{K} \mu_k^{m_k + \alpha_k -1} \\ &= \underbrace{\frac{B(m+\alpha)}{B(\alpha)}}_{\text{evidence }p(D|\alpha)} \,\underbrace{\mathrm{Dir}(\mu|m+\alpha)}_{\text{posterior }p(\mu|D,\alpha)} \end{align} $$
This equation is the equivalent of the Gaussian multiplication formula for discrete data. Note that the evidence is a scalar normalizer for given observations $m$ and pseudo-observations ("prior" observations) $\alpha$.
Next, we need a prior for the parameters $\mu = (\mu_1,\mu_2,\ldots,\mu_K)^T$.
In the binary coin toss example, we used a beta distribution that was conjugate with the binomial and forced us to choose prior pseudo-counts.
The generalization of the beta prior to $K$ parameters $\{\mu_k\}$ is the Dirichlet distribution:
$$p(\mu|\alpha) = \mathrm{Dir}(\mu|\alpha) = \frac{\Gamma\left(\sum_k \alpha_k\right)}{\Gamma(\alpha_1)\cdots \Gamma(\alpha_K)} \prod_{k=1}^K \mu_k^{\alpha_k-1} $$
where $\Gamma(\cdot)$ is the Gamma function.
The Gamma function can be interpreted as a generalization of the factorial function to the real ($\mathbb{R}$) numbers. If $n$ is a natural number ($1,2,3, \ldots $), then $\Gamma(n) = (n-1)!$, where $(n-1)! = (n-1)\cdot (n-2) \cdot 1$.
As before for the Beta distribution in the coin toss experiment, you can interpret $\alpha_k$ as the prior number of (pseudo-)observations that the die landed on the $k$-th face.
Consider a toss with a $K$-sided die. We use a one-hot coding scheme, i.e., the outcome is encoded as
$$x_{k} = \begin{cases} 1 & \text{if the throw landed on $k$-th face}\\ 0 & \text{otherwise} \end{cases} \,.$$
Assume the probabilities
$$p(x_{k}=1) = \mu_k \quad \text{with } \mu_k \geq 0 \text{ and }\sum_k \mu_k = 1 \,.$$
The data generating distribution for one-hot encoded outcome $x = (x_{1},x_{2},\ldots,x_{K})$ (and $\mu = (\mu_1,\mu_2,\dots,\mu_k)^T$) is then given by
$$p(x|\mu) = \mu_1^{x_1} \mu_2^{x_2} \cdots \mu_K^{x_K}=\prod_{k=1}^K \mu_k^{x_k} \tag{B-2.26}$$
This generalized Bernoulli distribution is called the categorical distribution.
Show that
(a) the categorial distribution is a special case of the multinomial for $N=1$.
(b) the Bernoulli is a special case of the categorial distribution for $K=2$.
(c) the binomial is a special case of the multinomial for $K=2$.
Show that Laplace's generalized rule of succession can be worked out to a prediction that is composed of a prior prediction and data-based correction term.
(a) The probability mass function of a multinomial distribution is
$$ p(D_m|\mu) =\frac{N!}{m_1! m_2!\ldots m_K!} \,\prod_k \mu_k^{m_k}$$
over the data frequencies $D_m=\{m_1,\ldots,m_K\}$ with constraints that $\sum_k \mu_k = 1$ and $\sum_k m_k=N$.
Setting $N=1$, we see that $p(D_m|\mu) \propto \prod_k \mu_k^{m_k}$ with $\sum_k m_k=1$, making the sample-space one-hot coded. This is the categorical distribution.
(b) When $K=2$, the constraint for the categorical distribution takes the form $m_1=1-m_2$ leading to
$$ p(D_m|\mu) \propto \mu_1^{m_1}(1-\mu_1)^{1-m_1}$$
which is associated with the Bernoulli distribution.
(c) Plugging $K=2$ into the multinomial distribution leads to $p(D_m|\mu) =\frac{N!}{m_1! m_2!}\mu_1^{m_1}\left(\mu_2^{m_2}\right)$ with the constraints $m_1+m_2=N$ and $\mu_1+\mu_2=1$. Then plugging the constraints back in we obtain
$$ p(D_m|\mu) = \frac{N!}{m_1! (N-m1)!}\mu_1^{m_1}\left(1-\mu_1\right)^{N-m_1}$$
which is the binomial distribution.
Derivations are in the lecture notes.
(a)
$$p(x_n|\mu) = \prod_k \mu_k^{x_{nk}} \quad \text{subject to} \quad \sum_k \mu_k = 1 \,.$$
$$p(D|\mu) = \sum_k m_k \log \mu_k$$
where $m_k = \sum_n x_{nk}$.
(b)
$$\hat \mu = \frac{m_k}{N}\,,$$
which is the sample proportion.
Simple Bayesian and maximum likelihood-based density estimation for discretely valued data sets
Mandatory
These lecture notes
Optional
Bishop PRML book (2006), pp. 67-70, 74-76, 93-94
When doing ML estimation, we must obey the constraint $\sum_k \mu_k = 1$, which can be accomplished by a Lagrange multiplier. The constrained log-likelihood with Lagrange multiplier is then
$$\tilde{\mathrm{L}}(\mu) = \sum_k m_k \log \mu_k + \lambda \cdot \big(1 - \sum_k \mu_k \big)$$
The method of Lagrange multipliers is a mathematical method for transforming a constrained optimization problem to an unconstrained optimization problem (see Bishop App.E). Unconstrained optimization problems can be solved by setting the derivative to zero.
The outcomes $x_n$ are encoded as
$$x_{nk} = \begin{cases} 1 & \text{if the $n$-th throw landed on $k$-th face}\\ 0 & \text{otherwise} \end{cases}$$
and the likelihood function for $\mu$ is now
$$p(D|\mu) = \prod_n \prod_k \mu_k^{x_{nk}} = \prod_k \mu_k^{\sum_n x_{nk}} = \prod_k \mu_k^{m_k} \tag{B-2.29}$$
where $m_k= \sum_n x_{nk}$ is the total number of occurrences that the outcome landed on face $k$. The vector $m = (m_1,m_2, \ldots, m_K)^T$ is known as the count vector. Note that $\sum_k m_k = N$.
This distribution depends on the observations only through the ''observed'' counts $\{m_k\}$. For given counts $\{m_k\}$, $p(D|\mu)$ can be interpreted as a likelihood function for $\mu$.
$$\begin{align*} p(&x_{\bullet,k}=1|D) = \frac{m_k + \alpha_k }{ N+ \sum_k \alpha_k} \\ &= \frac{m_k}{N+\sum_k \alpha_k} + \frac{\alpha_k}{N+\sum_k \alpha_k}\\ &= \frac{m_k}{N+\sum_k \alpha_k} \cdot \frac{N}{N} + \frac{\alpha_k}{N+\sum_k \alpha_k}\cdot \frac{\sum_k \alpha_k}{\sum_k\alpha_k} \\ &= \frac{N}{N+\sum_k \alpha_k} \cdot \frac{m_k}{N} + \frac{\sum_k \alpha_k}{N+\sum_k \alpha_k} \cdot \frac{\alpha_k}{\sum_k\alpha_k} \\ &= \frac{N}{N+\sum_k \alpha_k} \cdot \frac{m_k}{N} + \bigg( \frac{\sum_k \alpha_k}{N+\sum_k \alpha_k} + \underbrace{\frac{N}{N+\sum_k \alpha_k} - \frac{N}{N+\sum_k \alpha_k}}_{0}\bigg) \cdot \frac{\alpha_k}{\sum_k\alpha_k} \\ &= \frac{N}{N+\sum_k \alpha_k} \cdot \frac{m_k}{N} + \bigg( 1 - \frac{N}{N+\sum_k \alpha_k}\bigg) \cdot \frac{\alpha_k}{\sum_k\alpha_k} \\ &= \underbrace{\frac{\alpha_k}{\sum_k\alpha_k}}_{\text{prior prediction}} + \underbrace{\frac{N}{N+\sum_k \alpha_k} \cdot \underbrace{\left(\frac{m_k}{N} - \frac{\alpha_k}{\sum_k\alpha_k}\right)}_{\text{prediction error}}}_{\text{data-based correction}} \end{align*}$$
(If you know how to do it shorter and more elegantly, please post in Piazza.)
This decomposition is the natural consequence of doing Bayesian estimation, which always involves a prior-based prediction term and a likelihood-based (or data-based) correction term that can be interpreted as a (precision-weighted) prediction error.
We consider IID data $D = \{x_1,x_2,\ldots,x_N\}$ obtained from tossing a $K$-sided die. We use a binary selection variable
$$x_{nk} \equiv \begin{cases} 1 & \text{if $x_n$ lands on $k$-th face}\\ 0 & \text{otherwise} \end{cases}$$
with probabilities $p(x_{nk} = 1)=\mu_k$.
(a) Derive the log-likelihood $\log p(D|\mu)$.
(b) Derive the maximum likelihood estimate for $\mu$.
This is actually a generalization of the conjugate relation that we found for the binary coin toss:
$$\begin{align*} \underbrace{\text{beta}}_{\text{posterior}} &\propto \underbrace{\text{binomial}}_{\text{likelihood}} \cdot \underbrace{\text{beta}}_{\text{prior}} \end{align*}$$
In the above derivation, we noticed that the data generating distribution for $N$ die tosses with data outcomes $D=\{x_1,\ldots,x_N\}$ only depends on the counts $m_k$:
$$p(D|\mu) = \prod_n \underbrace{\prod_k \mu_k^{x_{nk}}}_{\text{categorical dist.}} = \prod_k \mu_k^{\sum_n x_{nk}} = \prod_k \mu_k^{m_k} \tag{B-2.29}$$
Now consider a $K$-sided coin (e.g., a six-faced die (pl.: dice)). How should we encode outcomes? Two natural options present themselves:
$$x \in \{1,2,\ldots,K\} \,.$$
E.g., for $K=6$, if the die lands on the 3rd face, then $x=3$.
This coding scheme is called label (or index) encoding.
$$x = (x_1,\ldots,x_K)^T $$
where $x_k$ are binary selection variables, given by
$$x_k = \begin{cases} 1 & \text{if die landed on $k$th face}\\ 0 & \text{otherwise} \end{cases}$$
For instance, for $K=6$, if the die lands on the $3$-rd face, then $x=(0,0,1,0,0,0)^T$.
This coding scheme is called a 1-of-K or one-hot coding scheme.
It turns out that the one-hot coding scheme is mathematically more convenient!
$$\begin{align*} \hat{\mu}_k &= \arg\max_{\mu_k} p(D|\mu) \\ &= \arg\max_{\mu_k} p(D|\mu) \cdot \underbrace{\left.\mathrm{Dir}(\mu|\alpha)\right|_{\alpha=(1,1,\ldots,1)}}_{\text{uniform distr.}} \\ &= \arg\max_{\mu_k} \left.p(\mu|D,\alpha)\right|_{\alpha=(1,1,\ldots,1)} \\ &= \arg\max_{\mu_k} \left.\mathrm{Dir}\left( \mu | m + \alpha \right)\right|_{\alpha=(1,1,\ldots,1)} \\ &= \frac{m_k}{\sum_k m_k} = \frac{m_k}{N} \end{align*}$$
where we used the fact that the maximum of the Dirichlet distribution $\mathrm{Dir}(\{\alpha_1,\ldots,\alpha_K\})$ is obtained at $(\alpha_k-1)/(\sum_k\alpha_k - K)$.
We can obtain the maximum likelihood estimate for $\mu_k$ based on $N$ throws of a $K$-sided die within the Bayesian framework by letting the prior for $\mu$ approach a uniform distribution. For a Dirichlet prior $\mathrm{Dir}(\mu | \alpha)$, this corresponds to setting $\alpha \rightarrow (1, 1, \dots, 1)$.
Prove for yourself that
$$\begin{align*} \hat{\mu}_k &= \arg\max_{\mu_k} p(D|\mu) = \frac{m_k}{N}\,. \end{align*}$$